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ABSTRACT 

This  paper  describes  a  diagnostic  objective  analysis  procedure  designed 
to  derive  dynamically  consistent  three-dimensional  fields  of  horizontal 
wind  velocity,  geopotential  and  temperature  over  limited-area  subsynoptic 
and  mesoscale  domains.  These  analyses  provide  data  for  initial  and  boundary 
conditions  for  the  fine-mesh  model  tests  over  the  United  States  and  Europe 
made  at  Fleet  Numerical  Weather  Central  (FNWC)  and  described  by  Anthes 
(1978).   Objective  analyses  are  produced  on  nine  standard  pressure  levels 
over  70  x  70  unstaggered  60-km  Cartesian  grids  on  Lambert  conformal  map  pro- 
jections covering  the  United  States  and  western  Europe  including  the 
western  Mediterranean  Sea. 

The  analysis  procedure  utilizes  first-guess  fields  of  geopotential, 
temperature,  relative  humidity,  geostrophic  wind  and  sea-level  pressure 
taken  from  operational  analyses  on  the  FNWC  northern  hemispheric  63  x  63 
381-km  grid.   Wind  observations  originating  from  rawinsondes,  aircraft, 
pilot  balloons  and  satellite-derived  cloud  motions  modify  the  geostrophic 
first-guess  wind  analyses  on  pressure  surfaces  through  a  successive 
correction  scheme  consisting  of  two  scans  using  ellipsoidal  (three-dimen- 
sional) weighting  functions.   Following  the  wind  analysis,  we  infer  geo- 
potential and  temperature  from  the  rotational  part  of  the  wind  with  a  non- 
linear form  of  the  balance  equation.   This  procedure  ensures  that  the  wind 
and  mass  fields  are  consistently  related  through  a  simplified  set  of  the 
equations  of  motion. 

We  present  detailed  results  from  two  cases  and  error  statistics  from 
thirty  cases  occurring  from  December  1977  through  April  1978.   Typical 
root-mean-square  (RMS)  differences  between  first-guess  and  balanced  analyses 
of  geopotential  and  temperature  are  ~20  m  and  2°C,  while  RMS  vector  differ- 
ences between  analyzed  and  balanced  winds  are  ~5  m  s"-*-,  a  reasonable  order 
of  magnitude  for  the  divergent  component  of  the  wind.   Discrepancies  between 
the  first-guess  and  balanced  temperature  fields,  largest  near  the  ground 
and  at  the  tropopause,  appear  to  be  a  consequence  of  ignoring  the  divergent 
part  of  the  wind  when  solving  for  the  geopotential  and  truncation  errors 
when  inferring  temperature  from  geopotential. 
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1.   Introduction 

The  development  and  implementation  of  operational  numerical  weather  pre- 
diction models  have  required  the  concomitant  application  of  objective  analysis 
schemes  to  provide  data  for  initial  conditions  and  verification,  and,  in  the 
case  of  limited-area  models  in  a  research  mode,  time-dependent  boundary  condi- 
tions.  The  initial  conditions  for  a  primitive  equation  model  not  only  should 
reproduce  the  features  in  the  flow  within  the  scales  of  motion  resolvable  by 
the  model,  but  should  imply  local  time  tendencies  of  the  dependent  variables 
characteristic  of  the  analyzed  features.   Failure  to  satisfy  the  latter  require- 
ment leads  to  gravity-wave  noise  early  in  the  forecast. 

Independent  analyses  of  the  wind  and  mass  (specified  by  the  surface  pres- 
sure and  three-dimensional  temperature  patterns)  invariably  lead  to  gravity- 
wave  generation  in  models  because  these  fields  are  mutually  related  through 
the  equations  of  motion.   Therefore,  however  carefully  one  prepares  indepen- 
dent analyses  of  the  wind  and  mass  field,  small  inconsistencies  between  these 
fields  will  result  in  spurious  initial  accelerations  typical  of  gravity  waves. 
The  traditional  technique  for  obtaining  mutually  consistent  initial  analyses 
for  mid-latitude  synoptic-scale  models  has  been  to  derive  the  wind  field  from 
the  mass  field  through  a  diagnostic  relationship  such  as  the  geostrophic  wind 
law  or  a  diagnostic  form  of  the  divergence  equation  known  as  the  balance 
equation. 

This  paper  explores  the  feasibility  of  utilizing  the  reverse  procedure  of 
inferring  the  mass  field  from  the  wind  field,  which  is  standard  procedure  in 
tropical  regions,  for  meso-  and  subsynoptic-scale  motions  in  mid-latitudes. 
Given  analyses  of  sea-level  pressure,  boundary  geopotential,  and  the  rota-      , 
tional  part  of  the  wind,  we  derive  interior  geopotential  and  temperature 


using  the  nonlinear  balance  equation.   The  scheme  was  tested  in  15  cases 
over  North  America  and  15  cases  over  Europe  and  North  Africa  from  December 
1977  through  April  1978.   Furthermore,  the  same  analysis  procedure  generated 
initial  and  boundary  data  for  the  version  of  The  Pennsylvania  State  Univer- 
sity mesoscale  model  (Anthes  and  Warner,  1978)  tested  at  the  Naval  Environ- 
mental Prediction  Research  Facility  (NEPRF)  and  Fleet  Numerical  Weather  Cen- 
tral (FNWC)  in  early  1978.   The  results  of  thirty-two  forecasts  over  the 
eastern  two-thirds  of  the  United  States  and  western  Europe,  including  the 
western  Mediterranean  Sea,  are  given  in  a  separate  report  (Anthes,  1978). 
2.   Analysis  procedure 

The  analysis  procedure  comprises  three  main  steps:  obtaining  first-guess 
fields  of  geopotential ,  temperature,  relative  humidity  and  horizontal  wind 
velocity;  producing  a  mesoscale  analysis  of  the  rotational  part  of  the  observed 
winds;  and  inferring  geopotential  and  temperature  from  the  nonlinear  balance 
equation. 

The  analyses  are  prepared  on  70  x  70  60-km  unstaggered  Cartesian  grids 
on  Lambert  conformal  map  projections  for  Europe  and  North  Africa  (grid  cen- 
tered at  45''N,  5°E)  and  North  America  (grid  centered  40°N,  90°W)  .   Figs,  la 
and  lb  show  the  geographical  extent  of  the  domains.   The  interior  region  of 
the  70  X  70  grids,  delineated  by  the  inner  box,  marks  the  50  x  50  60-km  model 
forecast  domain.   We  perform  the  analyses  on  the  larger  domain  to  minimize 
the  effect  of  boundary  wind  errors  on  the  model  initial  conditions  (Elsberry 
and  Ley,  1976).   This  rationale  is  justified  in  the  Appendix. 

a.   Determination  of  the  first-guess  fields 

We  obtain  first-guess  fields  of  sea-level  pressure,  geopotential ,  tem- 
perature, and  vapor  pressure  (which  is  immediately  converted  to  relative 
humidity)  on  a  rectangular  subsection  of  the  FNWC  63  x  63  381-km  northern 
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hemispheric  polar  stereographic  grid  that  includes  the  70  x  70  analysis  domain. 
Geopotential  and  temperature  are  taken  from  analyses  at  925,  850,  700,  500, 
400,  300,  250,  200,  150  and  100  mb;  vapor  pressure  is  obtained  from  a  12-h 
forecast  of  the  FNWC  operational  model  (Kesel  and  Winninghoff,  1972)  valid  at 
the  analysis  time  at  the  same  pressure  levels  at  and  below  300  mb.   After  com- 
puting geostrophic  winds  on  the  FNWC  subset  grid  at  all  pressure  levels,  we 
interpolate  all  the  fields  to  the  70  x  70  analysis  domain  using  cubic  splines 
under  tension  (Cline,  1974)  with  the  tension  factor  equal  to  one.   Although 
the  first-guess  moisture  analysis  does  not  contain  mesoscale  details,  we  leave 
it  unmodified  during  the  remainder  of  the  analysis  procedure  and  neglect  its 
effects  on  the  mass  and  wind  fields.   Nevertheless,  the  need  for  a  detailed, 
high-resolution  moisture  analysis  for  initializing  mesoscale  models  justifies 
additional  attention  to  this  problem.   Finally,  we  note  that  the  procedure 
for  obtaining  the  first-guess  fields  requires  about  100  s  of  execution  (CP) 
time  on  FNWC ' s  CDC  CYBER-175  computer, 
b.   Analysis  of  observed  winds 

After  obtaining  the  first-guess  analyses  on  the  70  x  70  analysis  grid, 
we  perform  a  three-dimensional  objective  analysis  of  observed  winds  on  pres- 
sure levels  using  a  procedure  developed  by  Barker  at  NEPRF  (personal  communi- 
cation) .   The  scheme  consists  of  modifying  the  geostrophic  first-guess  analy- 
ses by  performing  two  scans  with  ellipsoidal  (three-dimensional)  weighting 
functions,  which  allows  one  to  specify  the  minimum  scales  in  the  final  wind 
analysis  (~240  km  =  4As  in  this  case) . 

Wind  observations  originate  from  rawinsondes,  aircraft,  pilot  balloons, 
and  satellite-observed  cloud  motions.   The  wind  observations  are  not  corrected 
for  balloon  drift  or  observation  time  differing  slightly  from  the  analysis 
time.   Finally,  the  wind  analysis  scheme  is  coded  so  that  a  horizontal  vector 


wind  observation  is  contained  within  one  word  of  memory,  so  it  can  be  treated 
computationally  as  a  scalar.   For  ten  pressure  levels  on  the  70  x  70  domain, 
the  wind  analysis  requires  about  120  s  of  execution  time, 
c.   Solution  for  the  mass  field 

Given  detailed  wind  analyses  and  boundary  geopotential  on  the  70  x  70 
analysis  domain,  we  derive  the  interior  geopotential  and  temperature  fields 
with  the  nonlinear  balance  equation  and  the  hydrostatic  assumption. 

We  choose  to  deviate  from  standard  mid-latitude  practice  of  inferring 
wind  from  mass  (Shuman,  1957)  because  the  noise  in  the  observed  height  field 

(~10  m)  implies  an  unacceptable  noise  level  in  the  geostrophic  wind  of  16.3 

-1  -A  -1 

m  s   for  a  60-km  grid  spacing  and  the  Coriolis  parameter  f  =  10   s   .   The 

same  problem  occurs  in  tropical  analysis  of  synoptic-scale  motions  where  a 
decrease  in  f  compensates  for  the  increased  grid  spacing.   Therefore,  despite 
the  risk  of  incurring  errors  from  uncertainties  in  wind  observations, 
especially  in  situations  where  winds  are  strong  and  at  high  levels,  causing 
low  elevation  angles  of  pilot  balloons  or  rawinsondes,  we  favor  inferring 
mass  from  the  winds  (Fankhauser ,  1974;  Anthes,  1976). 

Considerations  of  a  mathematical  nature  provide  an  additional  argument 
for  starting  with  the  winds.   The  usual  procedure  of  solving  the  nonlinear 
balance  equation  (see  (12))  for  the  stream  function  from  geopotential  requires 
solution  of  an  elliptic  partial  differential  equation,  which  has  unique  solu- 
tions provided  that 

C^>-f.  (1) 

where   C   is  the  relative  geostrophic  vorticity  on  a  pressure  surface  (Haltiner, 
g 

1971) .   While  the  constraint  in  (1)  is  usually  not  critical  in  synoptic-scale 
analyses,  it  does  not  permit  realistic  mesoscale  resolution  of  anticyclonic 


shear  zones,  especially  those  associated  with  jet  streaks.   The  alternative 
of  deriving  mass  from  wind,  however,  only  requires  solution  of  a  Poisson 
equation,  which  is  always  elliptic  and  solvable  by  a  straightforward  procedure. 

The  first  stage  of  the  balancing  procedure  consists  of  deriving  the 
stream  function  from  the  rotational  part  of  the  observed  wind.   According  to 
Helmholtz's  theorem,  the  rotational  part  of  the  wind 

V^  =  k  X  mVi|^  ,  (2) 

where  ip   is  the  nondivergent  stream  function,   k  the  vertical  unit  vector, 
m  the  map-scale  factor,  and  V  =  i  -r 1-  j 


8x    -"  9y   ' 
P       P 
From  (2) ,  the  component  of  the  wind  normal  to  the  boundary  is 

V  =  -m  -r^  ,  (3) 

n      OS 


where  s   is  distance  along  the  boundary,  positive  when  taken  counterclockwise. 
Following  Phillips  (1958),  Brown  and  Neilon  (1961),  Bedient  and  Vederman  (1964), 
Anthes  (1976),  we  compute  the  boundary  stream  function,  4j  ,  by  adding  a 

D 

correction  factor  (constant  for  each  pressure  level)  that  forces  the  mean 

observed  velocity  component  normal  to  the  boundary,  v     ,  to  be  zero.   For- 

n 

mally,  we  write 

V  =  V     +  c  ,  (4) 

n    n 

where 


obs 


f  ^  ds  =  f  {-^ +  ^)  ds  =  0  , 

»/   m       *'   m       m 


(5) 


so  that 
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obs 

V 

n    J 
ds 
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The  magnitude  of  the  correction  factor,  c,  is  usually  less  than  1ms   for 
the  domains  we  consider.   From  (3),  (4)  and  (6),  we  conclude  that 


o-   f    - 

o     J        m 


iPgCs)  =  ipgCs^)  -  /   ^  ds  ,  (7) 


s 
o 

where   s    applies  at  the  southwest  corner  of  the  domain. 

Given  ip„(s)   from  (7),  we  determine  ^p      in  the  domain  interior  by 

o 

solving  the  Poisson  equation  implied  by  (2) 

V^i|;  =  \—  ,  (8) 

m 


where 

,      ^  ^    obs  ^   obs 

^obs    l.d      fV .    d  ,u . -,  .    . 

C    =  m  [-^  (- )  -  3-(- )J  .  (9) 

dx  m  dy  m 

Eq.  (8)  is  approximated  by  centered-in-space  differences  on  an  unstaggered  grid 
and  solved  with  an  efficient  direct  method  (Rosmond  and  Faulkner,  1976).   In 
finite-difference  form  the  Laplacian  operator  is 

V^*  .  (*^^^^  +  *,_^.  +  *.._^  +  ^,,^y  -  4*.,)/(As)^  ,    (10) 

where   i  and  j   are  indices  in  the  y   and  x  directions,  respectively.   One 
can  show,  however,  from  the  finite  difference  form  of  the  vorticity  (9)  that 

is  actually  the  consistent  form  ensuring  that  the  vorticity  of  the  balanced 
wind  equals  that  of  the  observed  wind.   We  choose  (10)  in  order  to  avoid  need- 
ing two  sets  of  boundary  conditions  while  noting  that  both  (10)  and  (11)  yield 
nondlvergent  winds. 
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The  second  stage  of  the  balancing  procedure  consists  of  solving  the 
nonlinear  balance  equation 

V^c^  =  V  .  am)   +   2m2j(|^  ,  |i)  (12) 

for  geopotential,  (^ .      In  (12),  J   is  the  Jacobian  operator,  and  we  include  the 
full  variation  of   f  ,  which  for  the  Lambert  conformal  projection  varies  in 
X   as  well  as   y  .   In  addition,  we  have  neglected  the  spatial  variation  of 
the  map-scale  factor.   Eq.  (12)  is  consistent  with  neglecting  the  horizontal 
velocity  divergence   D  ,  dD/dt  ,  and  friction  in  the  full  divergence  equation. 

The  right  side  of  (12)  is  known,  so  that  it  can  be  solved  by  the  same 
technique  used  for  (8)  with  the  first-guess  (^     fields  interpolated  to  the 
boundaries  of  the  70  x  70  grid  serving  as  boundary  conditions.   As  discussed  in 
the  Appendix,  the  lack  of  mesoscale  detail  in  (\>     on  the  boundaries  is  accept- 
able, since  for  a  given  amplitude  the  exponential  damping  of  the  influence  of 
a  boundary  harmonic  with  distance  into  the  domain  interior  is  faster  for  higher 
wavenumbers.   Accurately  specifying  the  large-scale  variation  of   (^   on  the 
boundaries  is  important,  and  we  believe  the  FNWC  analysis  is  sufficiently 
reliable  for  this  purpose. 

After  solving  (12)  ,  our  analysis  will  be  complete  once  the  tempera- 
ture is  determined.   Using  the  hydrostatic  assumption  on  the  derived  geo- 
potentials  and  first-guess  sea-level  pressure  fields,  we  determine  estimates 
of  temperatures  at  intermediate  levels  determined  by  the  logarithmic  mean 
pressure  between  the  analysis  levels.   Then  the  temperatures  are  interpolated 
linearly  with  respect  to  the  natural  logarithm  of  pressure  to  the  analysis 
levels.   The  balancing  procedure  requires  about  20  s  of  execution  time  compared 
to  a  total  of  240  s  for  the  entire  analysis  procedure. 
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3.   Examples  of  the  analysis  procedure 

In  order  to  illustrate  the  analysis  procedure,  we  examine  300-mb  winds  and 
500-mb  heights  and  temperatures  for  the  50  x  50  model  forecast  domains  over 
western  Europe  (Fig.  la)  and  the  United  States  (Fig.  lb)  for  0000  GMT  12  January 
1978  and  1200  GMT  25  January  1978,  respectively.   The  forecasts  based  on  these 
initial  conditions  are  discussed  by  Anthes  (1978).   Because  these  cases  feature 
large  height  and  temperature  gradients,  strong  winds  and  pronounced  streamline 
curvature,  they  are  atypical  of  situations  commonly  encountered,  and  thus  pro- 
vide a  stringent  test  of  the  assumptions  in  the  balancing  procedure, 
a.   European  case  of  0000  GMT  12  January  1978 

Fig.  2  shows  the  sea-level  pressure  pattern  over  western  Europe  that 
foreshadowed  cyclogenesis  in  the  Gulf  of  Genoa.   An  old  occluded  low  (992  mb) 
over  the  Netherlands  is  weakening  and  a  low  is  developing  in  the  Mediterranean 
Sea  south  of  France.   We  mention  that  the  sea-level  pressure  pattern  does  not 
exhibit  mesoscale  detail  because  it  is  the  unmodified  first-guess  interpolated 
from  the  FNWC  hemispheric  analysis. 

The  500-mb  temperatures  and  D  values  (h-5574m)  for  the  first-guess  and 
balanced  cases  for  0000  GMT  12  January  1978  are  displayed  in  Figs.  3a  and  3b, 
respectively.   The  temperatures  and  D  values  do  not  match  on  the  boundaries 
because  only  the  50  x  50  interior  points  of  the  70  x  70  analysis  domain  are 
shown.   Both  analyses  show  a  high-amplitude,  tilted  short-wave  trough  in  the 
height  fields  with  a  diffluent  region  over  the  area  of  incipient  cyclogenesis. 
The  balanced  analysis  contains  a  cutoff  low  above  the  occluding  surface  system, 
while  the  FNWC  first-guess  analysis  shows  a  sharp  trough. 

The  first-guess  temperature  analysis  (Fig.  3a)  contains  a  cold  dome 
(~  -35°C)  centered  over  the  Bay  of  Biscay  and  a  mid-tropospheric  frontal  zone 
at  the  base  of  the  trough  in  the  height  field.   The  balanced  analysis  (Fig.  3b) 
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Figure  2.   First-guess  sea-level  pressure  analysis  for  50  x  50 
European  domain  at  0000  GMT  12  January  1978.   Solid 
lines  are  deviations  of  sea-level  pressure  from 
1000  mb;  contour  interval  is  4  mb. 
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Figure  3.   Temperature  and  height  analyses  at  500  mb  for  50  x  50 

European  domain  at  0000  GMT  12  January  1978.   Solid  lines 
are  deviations  of  the  height  from  5574  m;  contour  inter- 
val is  60  m.   Dashed  lines  are  isotherms;  contour  inter- 
val is  5°C. 

a.  First-guess  analysis 

b.  Balanced  analysis 
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has  relaxed  the  temperature  gradient  in  the  frontal  zone  and  produced  a  small 
minimum  of  -40°C  over  the  west  coast  of  France,  which  is  about  -5°C  lower  than 
in  the  first-guess  analysis.   The  smoothing  of  the  frontal  zone  can  be  explained 
by  recalling  that  temperatures  are  computed  hydrostatically  over  layers  ~100  mb 
thick,  while  typical  mid-tropospheric  frontal  zones  are  50-100  mb  deep.   There- 
fore, the  region  of  high  static  stability  marking  the  frontal  zone  is  smoothed 
vertically,  and,  in  order  to  preserve  the  frontal  slope  and  thermal  wind 
balance,  horizontal  smoothing  results  as  well.   The  bias  towards  low  balanced 
temperatures  in  cyclonically  curved  systems  is  not  confined  to  this  case  and 
will  be  explained  later. 

Wind  analyses  at  300  mb ,  approximately  the  level  of  maximum  wind,  are 
shown  in  Figs.  4a,  4b  and  4c.   Fig.  4a  displays  the  result  of  interpolating 
geostrophic  winds  computed  on  the  FNWC  381-km  grid  to  the  60-km  analysis  grid. 
A  cyclonically  curved  jet  streak  appears  on  the  west  side  of  the  trough,  a 
broad  southwesterly  stream  of  air  covers  central  Europe,  and  a  weak  short-wave 
in  the  wind  pattern  is  located  over  Italy. 

The  observed  wind  analysis  (Fig.  4b) ,  while  reflecting  the  first-guess 
flow  pattern,  is  far  more  detailed.   In  particular,  the  70-m  s   jet  maximum 
over  the  northwestern  corner  of  Spain  is  10  m  s   stronger,  but  of  smaller-scale 
and  less  elongated.   Removing  the  observed  divergence  from  the  wind  (Fig.  4c) 
eliminates  much  of  the  small-scale  detail,  confirming  the  hypothesis  that  the 
magnitude  of  the  observed  divergent  part  of  the  wind  lies  within  the  noise 
level  of  the  observations  and  their  subsequent  analysis.   In  the  nondivergent 
analysis  (Fig.  4c),  the  jet  streak  is  weaker  than  in  the  observed  wind  analysis 
(65  m  s    compared  to  75  m  s   )  and  the  isotachs  are  more  elliptical.   Comparing 
the  nondivergent  analysis  (Fig.  4c)  with  the  first-guess  analysis  (Fig.  4a) 
shows  the  cross-jet  width  to  be  less  in  the  former.   In  particular,  we  note 
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Figure  4.   Wind  analyses  at  300  mb  for  50  x  50  European  domain 

at  0000  GMT  12  January  1978.   Solid  lines  are  isotachs; 
contour  interval  is  10  m  s~^.      On  wind  arrows,  flags, 
full  barbs  and  half  barbs  denote  speeds  of  50,  10  and 
5  m  s"-'-,  respectively. 

a.  First-guess  wind  analysis 

b.  Objective  analysis  of  observed  winds 

c.  Nondivergent  wind  analysis 
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that  the  area  covered  by  the  50  m  s    isotach  is  narrower.   In  addition  to  ren- 
dering a  more  realistic  jet  streak  than  that  in  the  first-guess  analysis,  the 
nondivergent  wind  analysis  preserves  the  short  wave  over  Italy. 

One  may  question  why  the  differences  in  wind  speeds  in  the  trough  between 
the  first-guess  and  balanced  analysis  are  not  greater  in  view  of  the  large 
cyclonic  streamline  curvature.   The  discrepancy  is  due  to  computing  the  first- 
guess  winds  over  a  2As  (~760  km)  interval  on  the  FNWC  grid,  where  truncation 
errors  weaken  the  supergradient  geostrophic  winds  to  the  extent  that  they  roughly 
equal  the  nondivergent  winds. 

Table  1  presents  detailed  statistical  summaries  of  the  winds,  heights 
and  temperatures  for  the  70  x  70  European  domain  at  0000  GMT  12  January  1978. 
Examination  of  the  root-mean-square  (RMS)  vector  differences  between  the  observed 
and  nondivergent  winds  shows  that  a  typical  value  of  3  m  s   is  on  the  order  of 
both  the  divergent  and  error  components  of  the  wind.   The  small  reduction  between 
the  RMS  vector  differences  and  the  RMS  vector  differences  corrected  for  the  mean 
(see  the  Appendix  for  the  definitions  of  the  statistical  quantities)  indicates 
that  the  bias  between  the  fields  is  small,  which  is  consistent  with  the  domain- 
averaged  divergent  (error)  component  being  negligible. 

The  statistics  for  the  first-guess  and  balanced  heights  exhibit  RMS 
differences  ranging  from  6  m  at  925  mb  to  22  m  at  100  mb.   When  normalized  by 
the  geometric  mean  of  the  standard  deviations  of  the  first-guess  and  balanced 
heights,  the  RMS  errors  are  seen  to  increase  in  proportion  to  the  increasing 
variance  in  the  height  fields  with  elevation.   The  magnitude  of  the  RMS  height 
differences  are  comparable  to  the  10-m  errors  Fankhauser  (1974)  quoted  for  the 
solution  of  the  complete  divergence  equation  for  geopotential  on  a  mesoscale 
grid  with  10-km  spacing.   The  temperature  statistics  indicate  that  the  discrep- 
ancies between  the  first-guess  and  balanced  analyses  are  largest  near  the  ground 
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Table  1.   Statistical  summaries  for  first-guess,  fg,  and  balanced,  b, 
winds,  V;  heights,  h;  and  temperatures,  T;  for  European 
domain  at  0000  GMT  12  January  1978.   Statistical  quantities 
are  defined  in  the  Appendix. 
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Table   1.       (concluded) 
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and  tropopause,  while  the  closest  agreement  is  in  the  middle  troposphere, 
b.   United  States  case  of  1200  GMT  25  January  1978 

Although  qualitative  differences  between  first-guess,  observed  and 
balanced  analyses  are  quite  similar  in  this  case  to  those  in  the  preceding 
European  case,  we  include  it  because  the  synoptic  situation  developed  into  the 
worst  blizzard  on  record  in  the  Ohio  Valley  with  a  central  pressure  of  958  mb 
at  Cleveland  at  1200  GMT  26  January  1978  and  hurricane-force  winds  near  the 
storm  center  (Wagner,  1978).   Anthes  (1978)  discusses  the  model  forecast  based 
on  the  initial  conditions  given  in  this  section. 

The  first-guess  sea-level  pressure  pattern  (Fig.  5)  exhibits  low  pres- 
sure centers  in  northern  Minnesota  and  the  Gulf  States.^  Comparison  with  the 
National  Meteorological  Center's  (NMC)  hand-analyzed  sea-level  pressure  chart 
(not  shown)  reveals  that  although  the  pressure  gradients  are  resolved  reason- 
ably well,  the  central  pressures  of  the  lows  (1000  mb  and  998  mb  for  the  nor- 
thern and  southern  systems,  respectively)  are  several  millibars  too  high  as 
expected  for  a  large-scale  analysis. 

The  first-guess  and  balanced  500-mb  charts  (Figs.  6a  and  6b)  show  a 
vigorous  short-wave  trough  over  Oklahoma  and  Texas  "phasing"  with  a  cutoff  low 
over  eastern  North  Dakota.   The  balanced  D  values  for  the  central  height  in 
this  low  are  about  60  m  higher  than  in  the  FNWC  first-guess  analysis.   The 
same  holds  for  the  Oklahoma-Texas  trough.   As  with  the  European  case,  the  cold 
bias  appears  in  the  balanced  analysis  in  regions  of  cyclonic  curvature,  while 
in  straight  flow  the  agreement  is  closer  as  revealed  by  the  relative  positions 
of  the  -20°C  isotherms  in  the  eastern  United  States  (compare  Figs.  6a  and  6b). 
Finally,  we  mention  that  the  500-mb  first-guess  temperatures  and  D-values  agree 
very  closely  with  the  NMC  analysis  and  observations  (not  shown) ,  so  for  veri- 
fication purposes,  we  will  consider  the  first-guess  height  and  temperature 
analyses  as  the  true  mass  field. 
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Figure  5.   First-guess  sea-level  pressure  analysis  for  50  x  50 
United  States  domain  at  1200  GMT  25  January  1978. 
Solid  lines  are  deviations  of  sea-level  pressure  from 
1000  mb;  contour  interval  is  4  mb. 
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Figure  6.   Temperature  and  height  analyses  at  500  mb  for  50  x  50 
United  States  domain  at  1200  GMT  25  January  1978. 
Solid  lines  are  deviations  of  the  height  from  5574  m; 
contour  interval  is  60  m.   Dashed  lipes  are  isotherms; 
contour  interval  is  5°C. 

a.  First-guess  analysis 

b.  Balanced  analysis 


27 


(6a) 


28 


I   I   I   I   I   I  1   1   I   I   I   I   I   I   I   I   I   I   I   [   I    I   1^1   I   I   I    1  ■» 


(6b) 


29 


The  300-mb  first-guess,  observed  and  nondivergent  wind  analyses  display 
characteristics  similar  to  their  European  counterparts.   Overall,  the  nondiver- 
gent wind  analysis  (Fig.  7c)  is  streakier  than  the  first-guess  patterns  (Fig. 
7a) ,  and  smoother  with  slightly  weaker  winds  than  the  observed  winds  (Fig.  7b) . 
As  in  the  European  case,  removing  the  divergence  from  the  observed  winds  effec- 
tively filters  small-scale  details.   Finally,  we  note  that  the  highly  asymmetric 
300-mb  wind  pattern  (Fig.  7c),  containing  jet  streaks  in  excess  of  65  m  s   to 
the  west  side  of  the  North  Dakota  closed  low  and  in  the  base  of  the  Texas 
trough,  is  conducive  to  subsequent  cyclogenesis  (Newton,  1958). 

In  order  to  provide  evidence  that  removing  the  divergence  from  the 
initial  wind  fields  and  deriving  a  consistent  analysis  of.    the  mass  field  are 
useful  in  reducing  low-frequency  oscillations  in  the  surface  pressure  field,  we 
present  a  time  series  of  surface  pressure  at  a  point  in  northwestern  Missouri 
for  forecasts  (Anthes,  1978)  beginning  with  1)  first-guess  wind,  geopotential 
and  temperature  (geostrophic  initialization)and  2)  nondivergent  wind,  balanced 
geopotential  and  temperature  (balanced  initialization)  at  1200  GMT  25  January 
1978.   The  results  in  Fig.  8  reveal  the  superiority  of  the  balanced  initial 
conditions  over  the  geostrophic  initial  conditions,  which  produce  unrealistic 
10-mb  amplitude  oscillations  during  the  first  8  h  of  the  forecast. 

Briefly  alluding  to  the  statistical  summary  of  this  case  in  Table  2, 
we  note  typical  RMS  vector  wind  differences  between  the  observed  and  nondivergent 
analyses  of  about  5ms   ,  RMS  differences  between  first-guess  and  balanced 
heights  no  more  than  20  m,  and  temperature  differences  on  the  order  of  2°C  except 
for  low  levels  (925  and  850  mb)  and  the  tropopause  (200  mb) .   As  we  will  show 
in  the  following  section,  the  biases  between  the  first-guess  and  balanced  mass 
fields  appearing  in  this  case  are  not  exceptional,  but  arise  in  many  additional 
cases  as  well. 
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Figure  7.   Wind  analyses  at  300  mb  for  50  x  50  United  States 
domain  at  1200  GMT  25  January  1978.   Solid  lines 
are  isotachs;  contour  interval  is  10  m  s"-*-.   On 
wind  arrows,  flags,  full  barbs  and  Half  barbs  denote 
speeds  of  50,  10  and  5  m  s"-*-,  respectively. 

a.  First-guess  wind  analysis 

b.  Objective  analysis  of  observed  winds 

c.  Nondivergent  wind  analysis 
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Table  2.   Statistical  summaries  for  first-guess,  fg,  and  balanced,  b, 
winds,  V;  heights,  h;  and  temperatures,  T;  for  United  States 
domain  at  1200  GMT  25  January  1978.   Statistical  quantities 
are  defined  in  the  Appendix. 
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Table    2.       (concluded) 
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4.   Discussion  of  errors  in  determining  the  mass  field 

We  tested  the  analysis  scheme  described  in  the  previous  sections  in  15 

European  and  15  United  States  cases  occurring  between  December  1977  and  April 

1978  (see  Table  3  for  times  and  dates)  in  order  to  reveal  systematic  errors. 

This  section  iS  devoted  to  analyzing  temperature  errors  in  order  to  evaluate 

the  ability  of  the  method  to  infer  temperature  patterns  from  the  rotational 

part  of  the  wind  field. 

Figures  9a  and  9b  contain  first-guess  and  balanced  temperature  soundings 

spatially  and  temporally  averaged  over  the  15  European  and  15  American  cases 

listed  in  Table  3.   Lapse  rates  of  temperature  in  both  sets  of  soundings 

approximately  conform  to  the  standard  atmosphere  (not  shpwn) .   The  tropopause, 

located  between  250  and  300  mb,  is  slightly  lower  than  standard  (225  mb)  as 

expected  in  winter. 

Both  sets  of  mean  soundings  exhibit  strikingly  similar  differences  between 

the  balanced,  T,  ,  and  first-guess,  T   ,  temperatures.   In  order  to  magnify  the 

time-space  averaged  differences,  we  plotted  them  versus  Inp   in  Fig.  10. 

Except  for  the  925-mb  level,  the  biases,  T,  -  T^  ,  are  nearly  identical.   The 

b    fg 

balanced  temperatures  are  warmer  in  the  upper  troposphere  and  lower  stratosphere 
(400-150  mb)  with  a  maximum  ~1°C  at  the  tropopause  (250  mb) ,  and  the  lower 
troposphere  (850  and  925  mb).   In  the  middle  troposphere  (700-400  mb) ,  the 
balanced  temperatures  are  lower  than  their  first-guess  counterparts  by  as  much 
as  0.5°C  at  500  mb,  which  is  consistent  with  the  results  cited  in  the  two 
individual  case  studies. 

Fig.  11  presents  the  RMS  differences  between  T   and  T    plotted  against 
ilnp  for  the  European  and  United  States  cases  (compare  with  Tables  1  and  2).   As 
with  the  bias  in  Fig.  10,  the  RMS  differences  are  similar  for  the  two  regions 
except  at  925  mb,  where  the  3.7°C  difference  for  the  American  cases  is  particularly 
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Table  3.   List  of  times  and  dates  of  European  and  United  States  cases. 
The  eight  digit  number  is  in  the  format  YYMMDDHH,  where  YY 
is  the  last  two  digits  of  the  year;  MM,  month;  DD ,  day,  and 
HH,  hour  (GMT). 


Europe 


United  States 


77121800 
77121812 
78011200 
78011A00 
78011500 
78011912 
78012000 
78012600 
78021312 
78021A00 
78030300 
78030312 
78040200 
78040212 
78040300 


77121812 

77121900 

77121912 

78012412 

78012500 

78012512 

78012600. 

78012612 

78012700 

78021300 

78021312 

78021400 

78021412 

78030300 

78030312 
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Figure  9.   Skew-T  vs.  £np  diagrams  showing  mean  first-guess  (solid 
lines)  and  mean  balanced  (dashed  lines)  soundings. 

a.  Domain-averaged  for  all  Mediterranean  (European)  cases 

b.  Domain-averaged  for  all  United  States  cases 
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Figure  10.   Differences  (°C)  between  mean  balanced  and  first- 
guess  soundings  shown  in  Fig.  9  vs.  ilnp.   Solid  line 
depicts  United  States  cases;  dashed  line  depicts 
Mediterranean  (European)  cases. 
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Figure  11.   Root-mean  square  differences  (°C)  between  balanced  and 
first-guess  temperatures  vs.  £np  for  all  United  States 
cases  (solid  line)  and  all  Mediterranean  (European) 
cases  (dashed  line) . 
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disconcerting.   In  addition  to  the  lowest  level,  errors  are  also  large  (~2.5°C) 
at  the  tropopause  and  least  (~1.2°C)  in  the  middle  troposphere.   The  magnitude 
of  these  RMS  differences  may  be  compared  against  typical  rawinsonde  temperature 
errors  of  around  1°C. 

The  results  of  Figs.  10  and  11  reveal  four  problems  requiring  attention: 
1)  the  positive  bias  at  the  tropopause;  2)  the  negative  bias  in  the  middle  tropo- 
sphere; 3)  the  positive  bias  at  925  mb;  and  4)  the  difference  between  the 
European  and  American  biases  at  925  mb. 

We  begin  addressing  these  issues  by  examining  individual  soundings  for  the 
United  States  case  of  1200  GMT  25  January  1978.   Fig.  12a  and  Table  4a  contain 
the  sounding  for  41.9°N,  78.4°W,  a  location  on  the  western  Pennsylvania-New  York 
border;  Fig.  12b  and  Table  4b  show  the  sounding  for  47.5°N,  102. 0°W  in  western 
North  Dakota.   In  addition  to  exhibiting  the  first-guess  and  balanced  soundings, 
Figs.  12a  and  12b  contain  soundings  computed  from  the  first-guess  height  field. 
The  "x's"  indicate  the  temperature  estimates  at  intermediate  levels  implied  by 
the  first-guess  height  field  through  the  hydrostatic  approximation,  which  are 
used  to  derive  the  temperature  curve  by  linearly  interpolating  back  to  mandatory 
levels . 

The  tropopause  in  the  first-guess  sounding  in  Fig.  12a  is  located  around 
200  mb,  and,  consistent  with  the  time-space  averaged  statistics,  colder  than 
the  balanced  sounding  or  that  computed  from  the  first-guess  height  field.   The 
warm  discrepancy  between  the  balanced  and  first-guess  tropopause  temperatures 
may  be  due  partly  to  the  interpolation  from  the  intermediate  to  standard  pres- 
sure levels,  which  is  unable  to  resolve  the  minimum  at  the  tropopause.   In  the 
case  of  an  analytically  specified  sounding  with  a  constant  lapse  rate  in  the 
troposphere  (y  =  0.6  Y,),  and  an  isothermal  lapse  rate  above  200  mb,  differences 
between  the  hydrostatically  computed  and  specified  soundings  were  less  than  0.1°C 
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Figure  12.   Skew-T  vs.  inp   diagrams  for  two  locations  at  1200  GMT 
25  January  1978.   Solid  line  depicts  sounding  computed 
from  first-guess  heights;  "x's"  denote  layer-averaged, 
hydrostatically  derived  temperature^  used  to  construct 
this  sounding.   Long  dashes  denote  first-guess  sounding; 
short  dashes  denote  balanced  sounding. 

a.  Sounding  at  41.9°N,  78.4°W 

b.  Sounding  at  47.5°N,  102. 0°W 
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(12a) 


47 


■100 


(12b) 
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Table  4.   Sounding  data  for  first-guess,  fg,  and  balanced,  b,  tempera- 
tures, T,  and  heights,  h,  for  1200  GMT  25  January  1978. 

a.  Sounding  at  41.9°N,  78.4°W. 

b.  Sounding  at  47.5°N,  102. 0°W. 


p 
(mb) 

fg 

^b 
(°C) 

b   fg 

fg 

(m) 

b   fg 

a.   925 

1.7 

2.1 

0.4 

706 

702 

-4 

850 

1.8 

0.6 

-1.2 

1384 

1380 

-4 

700 

-4.1 

-4.2 

-0.1 

2942 

2934 

-8 

500 

-19.5 

-20.6 

-1.1 

5530 

5513 

-17 

400 

-31.1 

-32.7 

-1.6 

7153 

-  7126 

-27 

..   300 

-46.2 

-47.3 

-1.1 

9123 

9086 

-37 

250 

-56.1 

-54.0 

2.1 

10303 

10270 

-33 

200 

-63.4 

-57.7 

5.7 

11694 

11682 

-12 

150 

-57.1 

-59.5 

-2.3 

13493 

13492 

-1 

100 

-59.7 

16025 

b.   925 

-25.9 

-11.7 

14.2 

742 

788 

46 

850 

-26.9 

-24.3 

2.6 

1350 

1407 

57 

700 

-29.8 

-30.7 

-0.9 

2751 

2810 

59 

500 

-34.6 

-40.7 

-6.1 

5109 

5129 

20 

400 

-39.9 

-41.6 

-1.7 

6651 

6637 

-14 

300 

-43.9 

-38.5 

5.4 

8623 

8597 

-26 

250 

-44.6 

-42.0 

2.6 

9872 

9856 

-16 

200 

-45.8 

-47.7 

-1.9 

11366 

11331 

-35 

150 

-45.4 

-48.0 

-2.6 

13273 

13228 

-45 

100 

-52.3 

15936 
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in  the  troposphere  (at  and  below  250  mb)  and  2.3°C  at  200  mb,  the  assumed  tropo- 
pause  level.   In  the  real  atmosphere  we  may  expect  this  difference  to  be  smaller 
since  the  change  in  lapse  rate  at  the  tropopause  is  usually  less  abrupt.   A 
related  source  of  error  at  the  tropopause  is  revealed  in  Fig.  12b,  where  the 
balanced  and  hydrostatically  computed  temperatures  are  quite  similar  but  warmer 
than  the  first-guess  temperature  at  300  and  250  mb.   This  latter  problem  suggests 
an  inconsistency  between  our  method  of  hydrostatically  computing  temperatures 
and  FNWC's  method  of  relating  heights  and  temperatures.   Whether  or  not  this 
.latter  problem  is  systematic  is  unknown. 

The  former  sources  of  error  are  a  result  of  vertical  truncation  and  would 
occur  even  if  the  balanced  heights  were  perfect.   A  final  line  of  speculation 
concerning  the  positive  bias  in  the  high  troposphere  can  be  attributed  to  errors 
in  the  balanced  heights.   Ignoring  vertical  truncation  errors,  we  can  relate 
the  discrepancy  in  the  temperature  to  the  discrepancy  in  the  thickness  hydro- 
statically so  that 

8h'    .RT  .  (13) 


9  £np 


where  the  prime  denotes  the  difference  between  a  balanced  and  first-guess  quan- 
tity [(   )'=(    ),  ~(    )^]«   According  to  (13),  temperature  differences 
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T'   are  positive  (negative)  above  (below)  minima  in  h'  .   For  the  sounding  in 
Fig.  12a,  h'   has  a  minimum  of  -37  m  at  300  mb  increasing  to  -1  m  at  150  mb 
(Table  4a),  which,  from  (13),  implies  T*  =  1.8°C.   On  the  other  hand,   h'   in- 
creases in  the  lower  troposphere  to  -4  m  at  850  mb  from  the  300  mb  value, 
implying  T*  =  -1.1°C  between  these  two  levels.   The  trend  is  similar  in  Fig. 
12b  and  Table  4b  where  balanced  heights  are  higher  in  the  lower  troposphere  and 
lower  in  the  high  troposphere  than  the  first-guess  values  (h'  =  59  m  at  700  mb, 
-26  m  at  300  mb  implying  T'  =  -3.4°C),  a  situation  that  can  explain  the  negative 
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bias  in  the  middle  troposphere  and  the  positive  bias  at  the  surface  and  tropo- 
pause.   This  trend  is  evident  in  the  domain- averaged  height  differences  for 
this  time  as  well  (see  Table  2). 

The  temperature  calculation  at  925  mb  is  extremely  sensitive  to  height 
errors  since  the  temperature  is  computed  over  a  thin  layer.   For  example,  for 
h'  =  10  m  and  the  sea-level  pressure  equal  to  1013  mb,  T'  =  3.8°C  (see  Table  2 
for  comparison  with  925  mb  height  and  temperature  statistics).   Moreover,  these 
temperature  errors  become  more  serious  in  the  vicinity  of  low  pressure  systems. 
If  the  sea-level  pressure  were  1000  mb,  h'  =  10  m,  then  T'  would  equal  4.4°C. 

The  order  of  magnitude  and  sign  of  the  46-m  height  discrepancy  at  925  mb 
for  the  North  Dakota  sounding  (Table  4b)  can  be  explained  by  the  divergence  and 
friction  terms  neglected  in  forming  the  balance  equation  (12).   If  the 
tilting  terms  are  unimportant,  the  bias  in  geopotential  (()'   is  approximately 
given  by 

V^4)'  =  D^  +4?  +  kD  ,  (14) 

dt 

where  D  is  the  horizontal  divergence  and  k  a  linear  friction  coefficient 

given  by   C  |v|/H  ,  where   C   is  a  drag  coefficient,   |v|   surface  wind  speed 

and  H  boundary  layer  depth.   A  crude  scale  analysis  allows  us  to  revise  (14) 

into 

(1)    (2)   (3) 

Here  L   is  the  wavelength  of  a  surface  cyclone,  and  V  a  phase  propagation 
speed  so  that  the  time  scale,  T  =  ^  ,  for  the  parcel  derivative  is  assumed 
to  be  due  to  the  local  variation  at  a  point.   For  L  ~  2000  km  (Fig.  5), 
V  ~  10  m  s"^,  k  ~  IxlO"^  s"-^  (Cjj  =  5x10"^,  H  =  500  m),and  D  ~  -IxlO"^  s" 
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(divergence  appropriate  for  medium-intensity  synoptic  systems  [Petterssen,  1956]), 
T  is  about  14  h  and 

(1)      (2)    (3) 
■  h'  ~  -2.5  m+5m+25m=  27.5  m  .  (16) 

From  the  above  result,  we  see  that  ignoring  frictional  effects  in  the  planetary 
boundary  layer  (PBL)  dominates  as  a  source  of  error.   We  can  use  the  results  of 
the  scale  analysis  to  explain  the  differences  in  the  925-mb  temperature  biases 
between  the  European  and  American  cases.   Although  we  do  not  offer  quantitative 
evidence,  inspection  of  the  individual  cases  reveals  that  the  American  cases 
are  dominated  by  deeper,  more  intense  low  pressure  systems,  where  ignoring  sur- 
face friction  in  the  presence  of  low  sea-level  pressure  magnifies  the  positive 
bias  between  the  balanced  and  first-guess  925-mb  temperatures.   Furthermore, 
we  note  that  the  925-mb  temperature  bias  is  negligible  for  the  United  States 
sounding  (Fig.  12a,  Table  4a)  over  a  flat  sea-level  pressure  field  (Fig.  5)  com- 
pared to  the  sounding  (Fig.  12b,  Table  4b)  where  the  cyclonic  curvature  is  strong. 

A  similar  calculation  may  be  performed  for  the  tropopause  level  (~300  mb). 
Taking  L  ~  4000  km  (Fig.  6a),  V  ~  40  m  s    (here  V  is  a  wind  velocity  so  that 
T   is  an  advective  time  scale),  k  =  0,  and  D  ~  1x10    s   ,  T   is  about   7  h  and 

(1)      (2)     (3) 
h'  ~  -10  m+-40m+0m=-50m.  (17) 

The  major  omission  in  this  case  is  the  time  derivative  of  the  divergence  follow- 
ing parcel  trajectories. 
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In  summary,  the  discrepancies  between  the  balanced  and  first-guess  tempera- 
ture analyses  can  be  attributed   to  1)  vertical  truncation  or  interpolation 
error  in  attempting  to  resolve  the  discontinuity  in  lapse  rate  at  the  tropo- 
pause,  2)  ignoring  the  divergent  part  of  the  wind  and  surface  friction  in  deter- 
mining geopotential,  and  3)  inconsistencies  between  our  method  and  FNWC's 
method  for  hydrostatically  relating  heights  and  temperatures. 
5.   Summary  and  suggestions  for  future  investigation 

Extensive  tests  of  solving  for  the  mass  field  from  the  rotational  part  of 
the  wind  field  over  a  large  number  (30)  of  cases  indicate  that  RMS  heights  and 
temperatures  can  be  resolved  within  20  m  and  2°C,  respectively.   Errors  in 
boundary  layer  temperature  may  be  larger  locally  on  account  of  neglecting 
friction  and  because  of  the  shallow  depth  of  the  layer.   The  above  figures  apply 
to  rectangular  mid-latitude  domains  on  the  order  of  4000  km  on  a  side  with  60-km 
grid  resolution.   The  boundary  geopotential  are  assumed  known  from  a  reliable 
large-scale  analysis. 

There  is  strong  evidence  that  including  an  analysis  of  the  divergent  part 
of  the  wind  and  surface  friction  would  result  in  an  improvement  in  determining 
the  temperature  field.   A  possibility  for  including  the  divergent  part  of  the 
wind  in  the  analysis  in  a  manner  consistent  with  filtering  gravity  waves  lies 
in  experimenting  with  the  vorticity  form  of  the  balance  equation, 

V'V 
\l^<p   =  C(C+f)  +  k»V(C+f)  X  V  -  V   (-=Y^)  +  V'F^   ,  (18) 

where  F   is  a  friction  term  and  we  have  assumed  that  ^r-  =  0  .   Saha  and 
~r  "'- 

Suryanarayana  (1971)  compared  (18)  with  F  =  0  with  the  nonlinear  balance 


This  list  excludes  the  effect  of  errors  in  the  rotational  part  of  the  wind 
on  the  derived  geopotential.   We  have  assumed  these  errors  are  random  and  do 
not  contribute  to  the  bias  in  the  derived  temperatures. 
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equation  (12) ,  the  linear  balance  equation  and  quasi-geostrophic  form  of  the 
balance  equation,  when  solving  for  geopotentlal  in  the  tropics.   They  found 
(18)  to  give  results  slightly  superior  to  the  nonlinear  and  linear  forms, 
while  the  quasi-geostrophic  version  yielded  results  inferior  to  the  three 
former  balance  equations. 

Despite  the  improvement  in  the  temperature  patterns  expected  from  includ- 
ing the  divergent  part  of  the  wind  and  surface  friction  in  the  analysis, 
determining      the  former  is  a  vexing  problem  because  it  is  masked  by  the 
errors  in  the  wind  data.   The  method  devised  by  Endlich  (1967) ,  which  adjusts 
a  vector  field  to  imposed  values  of  vorticity  or  divergence  may  offer  a  viable 
alternative.   Adjustment  of  the  horizontal  wind  components  to  satisfy  boundary 
conditions  on  the  vertical  velocity  (O'Brien,  1970;  Fankhauser,  1974)  is  an 
additional  possibility. 

Finally,  we  point  out  that  the  gain  in  accuracy  with  our  approach  of 
solving  for  mass  from  winds  results  from  integrating  the  velocity  to  obtain 
stream  function  and  geopotentlal,  a  process  that  tends  to  smooth  small-scale 
observational  errors.   We  take  a  loss  in  accuracy,  however,  when  we  vertically 
differentiate  the  geopotentlal  to  obtain  temperature,  because  differentiating 
real  data  tends  to  introduce  noise.   Despite  such  inherent  limitations,  our 
method  for  inferring  the  mass  field  is  competitive  with  existing  analysis 
schemes,  and  has  the  potential  for  further  improvement. 
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APPENDIX 

1.   Comments  on  solutions  to  Poisson  equations 

When  solving  for  the  nondlvergent  wind  and  the  balanced  geopotentlal, 
we  must  solve  a  Poisson  equation 

V^(|.(x,y)  =  F(x,y)  (19) 

given  (p     on  the  lateral  boundaries.   We  proceed  to  Interpret  the  analytic 
solution  to  (19)  In  order  to  help  understand  the  numerical  results  obtained 
through  the  analysis  procedure. 

For  a  rectangular  domain  with  0  <^  x  j<  L  and  0  _<  y  £  D,  the  solution 
to  (1)  is 

D 

o  ^-,    slnh  :r—   (L-x)  „  r  ^    1 

^(-'y)=|E    hi --"^W       ^(0,y')sln^     dy') 

n=l   "^"^  "F"  J 


D 

<»        ,   .    nir 


?   -v        slnh  —  X  TT  r 

bZL       .        nTTL  "^"VW       ^a,y')sln^ 

n=l   ^^^^  -T-  J 


^■"^     dy') 


00  nIT 

slnh  -—  (D-y)  r  , 

.    ,   nTTD "^^"irW       ^(^'.0)sxn—     dx') 


n=l  L 

o 


oo  nTTy 


slnh      ^  ,  , 

smh 


/ 


n=l  L 

o 

D  L 
nTTx      .      mTTy 


'^i    i    Z  )  W  )      (//F(=<'.y)sinBM'3.„^'..M.., 


4LD 

n=l  m=l  o   o 

(20) 
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The  first  four  terms  on  the  right  are  the  contributions  of  each  lateral 
boundary  condition  to  the  solution.  The  final  term  is  the  component  of 
the  solution  due  to  the  forcing  in  the  interior. 

Because  (19)  is  linear,  we  may  interpret  (20)  for  a  single  harmonic 
without  loss  of  generality.   We  see  that  the  influence  of  the  right  lateral 
boundary  condition,  c{)(L,y)  [the  second  term  on  the  right  of  (20)]  damps 
exponentially  with  distance  from  the  boundary  at  a  rate  determined  by  the 
wavenumber  of  the  harmonic  representing  (f)(L,y).   This  point  is  illustrated 
in  Fig.  13,  where  the  term  controlling  the  damping  of  the  boundary  condi- 
tion 4)(L,y)  is  graphed  as  a  function  of  wavenumber  — —  and  normalized 
distance  —  .   For  a  given  distance  from  the  boundary,  for  increasing  wave- 
number,  the  value  of  the  function  decreases. 

The  implications  of  the  preceding  arguments  are: 

1)  The  influence  of  the  boundary  conditions  on  the  solution  decreases 
exponentially  with  distance  from  the  boundaries. 

2)  It  is  important  to  specify  the  large-scale  (low  wavenumber)  varia- 
tion of  the  boundary  conditions  accurately,  since  amplitude  errors  for  low 
wavenumber  components  of  the  boundary  conditions  damp  less  rapidly. 

The  preceding  analysis  is  modified  slightly  when  a  discrete  domain  is 
considered  (R.  T.  Williams,  personal  communication).   The  Laplace  equation 
that  is  solved  in  the  continuous  case  to  isolate  the  influence  of  the 
boundary  conditions  becomes 

^j+lk    jk+1    jk-1    J-lk     jk 

In  (21) ,  which  applies  to  a  grid  with  uniform  spacing  such  that  x  =  jAs 
and  y  =  kAs,  one  can  show  that 

^  yjAs  iAkAs  ,__, 

d)  ,  =  A  e     e  (22) 
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the  boundary  condition  at  x/L  =  1  as  a  function  of 
distance  normal  to  the  boundary  and  wavenumber.   See 
text  for  definition  of  symbols. 
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is  a  solution  to  (21)  provided  that 


cosh  yAs  =  2  -  cos  AAs  .  (23) 

271 

In  (22),  A  is  a  complex  constant,  X  =  - —  ,  where  L    is  wavelength  in 

L  y 

y        27T 

the  y  direction  such  that  L  =  KAs  ,  and  \i  =  - —  with  L   wavelength  in 

y  L       X        * 

X 

the  X  direction. 
We  can  write 


L  =  rL  (24) 

X     y 

or  M  =  -  (25) 


Fig.  14  is  a  graph  of  r  versus  K  as  determined  by  (23).   For  short 
wavelengths  with  respect  to  the  grid  the  damping  is  slower  than  in  the 
continuous  case  with  the  largest  difference  for  2As  waves  (r  =  1.8).   For 
waves  longer  than  lOAs  the  differences  between  the  continuous  and  discrete 
cases  are  negligible.   The  conclusions  drawn  from  the  continuous  case  apply 
to  the  discrete  case  provided  that  the  correction  for  short  waves  shown 
in  Fig.  14  is  considered. 
2.   Definitions  of  statistical  quantities 

For  a  quantity  a  we  denote  the  mean  over  the  domain  as  a  .   The 
mean-square  difference  (MSD)  between  two  sets  of  a  data  is 


MSD(a^  2^  =  Tctj^-Ct^T^  , 


(26) 


where  1  and  2  identify  their  respective  data  sets.   The  MSD  corrected  for 
the  mean  is  denoted  as 


MSDC(a^  ^)   =  I (a^  -  a^)  -  (a^  -  a^) | 2  (27) 
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and  the  variance  of  a  is  written  as 


s^  =  |a  -  a|2  .  (28) 


The  above  definitions  apply  to  scalar  or  vector  quantities  provided  that 
one  invokes  the  appropriate  definition  of  the  absolute  value  operator. 
If  one  defines  a  dimensionless  skill  score 

MSD(a,  _) 


S^  =  1^21  ^  (29) 


^1  "2 


it  can  be  shown  that 

2    ^    ^2         1^1  "  ^2'^ 

S^  =  —  +  —  -  2r  +  -^ ^  ,  (30) 

^2   "1  "1  ^2 


where   r   is  the  correlation  between  a   and  a    ,  so  that 


(a^  -  a^)  •  (a^  -  a^) 
r  = 


^  "2 


Defining  a  skill  score  corrected  for  the  mean  as 


one  may  show  that 


so  that 


(31) 


-   MSDC  (a^  _) 

S  ^  =  ^-^^  ,  (32) 

^1  ^2 


S  2  ^  i  +  ^  .  2r  (33) 

c    S2   s^ 


s2  ^  S  2  _^_1 2_  ^  (34) 

s   s 
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In  order  to  interpret  the  skill  score  S,  defined  in  (29),  we  let  the 
subscripts  1  and  2  denote  the  predicted  (or  derived)  field  and  the  verifi- 
cation field,  respectively.   For  a  given  MSD,  dividing  by  the  product  of 
the  standard  deviations  of  both  fields  attributes  greater  skill  (smaller 
value  of  S)  when  the  variation  of  the  fields  about  their  mean  is  larger, 
and  reproducing  the  verification  field  is  inherently  more  difficult. 
According  to  (30),  S  is  equal  to  zero  (a  perfect  score)  when  the  standard 
deviations  of  the  derived  and  verification  fields  are  equal,  their  means 
are  equal,  and  the  correlation  between  them  is  equal  to  one.   Moreover, 
we  observe  that  the  sum  of  the  first  two  terms  in  (30)  is  a  minimum  (equal 
to  2)  when  s^  =  s„  .   Therefore,  (29)  does  not  reward  a  conservative 
estimate  of  the  verification  field  (s^  <  s„) ,  while  the  MSD  itself  tends 
to  be  smaller  for  a  smoother  predicted  pattern. 


I  ^-,  r 
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